home *** CD-ROM | disk | FTP | other *** search
/ Shareware Grab Bag / Shareware Grab Bag.iso / 007 / c_math.arc / TAN.C < prev    next >
Encoding:
C/C++ Source or Header  |  1983-07-02  |  1.4 KB  |  74 lines
  1. /*
  2.     floating point tangent
  3.  
  4.     A series is used after range reduction.
  5.     Coefficients are #4285 from Hart & Cheney. (19.74D)
  6. */
  7.  
  8. #include <errno.h>
  9. #include <math.h>
  10.  
  11. int    errno;
  12. static double invpi      = 1.27323954473516268;
  13. static double p0     = -0.1306820264754825668269611177e+5;
  14. static double p1      = 0.1055970901714953193602353981e+4;
  15. static double p2     = -0.1550685653483266376941705728e+2;
  16. static double p3      = 0.3422554387241003435328470489e-1;
  17. static double p4      = 0.3386638642677172096076369e-4;
  18. static double q0     = -0.1663895238947119001851464661e+5;
  19. static double q1      = 0.4765751362916483698926655581e+4;
  20. static double q2     = -0.1555033164031709966900124574e+3;
  21.  
  22. double
  23. tan(arg)
  24. double arg;
  25. {
  26.     double modf();
  27.     double sign, temp, e, x, xsq;
  28.     int flag, i;
  29.  
  30.     flag = 0;
  31.     sign = 1.;
  32.     if(arg < 0.){
  33.         arg = -arg;
  34.         sign = -1.;
  35.     }
  36.     arg = arg*invpi;   /*overflow?*/
  37.     x = modf(arg,&e);
  38.     i = e;
  39.     switch(i%4) {
  40.     case 1:
  41.         x = 1. - x;
  42.         flag = 1;
  43.         break;
  44.  
  45.     case 2:
  46.         sign = - sign;
  47.         flag = 1;
  48.         break;
  49.  
  50.     case 3:
  51.         x = 1. - x;
  52.         sign = - sign;
  53.         break;
  54.  
  55.     case 0:
  56.         break;
  57.     }
  58.  
  59.     xsq = x*x;
  60.     temp = ((((p4*xsq+p3)*xsq+p2)*xsq+p1)*xsq+p0)*x;
  61.     temp = temp/(((1.0*xsq+q2)*xsq+q1)*xsq+q0);
  62.  
  63.     if(flag == 1) {
  64.         if(temp == 0.) {
  65.             errno = ERANGE;
  66.             if (sign>0)
  67.                 return(HUGE);
  68.             return(-HUGE);
  69.         }
  70.         temp = 1./temp;
  71.     }
  72.     return(sign*temp);
  73. }
  74.